Vamos a empezar explicando datos espaciales de Sistemas de Información Geográfica (SIG, GIS), para:

Información geográfica (codificada en un Sistema de computación, en este caso R)

Tipos de datos espaciales

setwd("~/Desktop/VICENÇ MUT inicio/OHW Intercoonecta/OHW 2025/Taller intermedio Oct 2025/Marina/Datos espaciales de métricas climáticas")

Fuentes de datos on-line

EmodNET → https://emodnet.ec.europa.eu/en/bathymetry Natural Earth → http://www.naturalearthdata.com USGS → https://earthexplorer.usgs.gov NASA → https://sedac.ciesin.columbia.edu Open Topography → https://opentopography.org UNEP → http://geodata.grid.unep.ch FAO → http://www.fao.org/geonetwork/srv/en/main.home

Ahora pasamos a dibujar/plotear y analizar datos espaciales (y temporales) mediante variables climáticas como la temperatura superficial del mar

Aquí se proporciona el código de R necesario para calcular varias métricas climáticas del paquete de R, VoCC, desarrollado en el artículo García Molinos et al. (2019) y aplicado en el artículo Sanz-Martín et al. (2024), además de en otros muchos artículos. Consúltalos para obtener información más detallada.

García Molinos et al. (2019): https://besjournals.onlinelibrary.wiley.com/doi/epdf/10.1111/2041-210X.13295

Paquetes importantes a usar para trabajar SIG en R raster: importar, exportar y analizar datos ráster sp: analizar datos espaciales (vectoriales)

Para este tutorial, necesitamos los siguientes paquetes que debemos conseguir:

También necesitamos un conjunto de datos de temperatura superficial que se puedan descargar de servidores de datos satelitales. En este caso utilizaremos Copernicus (página web https://data.marine.copernicus.eu/product/SST_MED_SST_L4_REP_OBSERVATIONS_010_021/description) pero que también están disponibles en el repositorio de datos de Intercoonecta Github (https://github.com/xxx completar cuando esté subido).

Cálculo de promedios en función de la frecuencia temporal de nuestros datos

Agregación de las celdas

Cálculo de la tendencia temporal

# para calcular el promedio de todos los años. Pero eso no lo utilizaremos más adelante. 
mean_all <- mean(yr_sst)
plot(mean_all, main ="Promedio de SST (ºC) 2000 a 2020") 

# para calcular la tendencia climática local 
?tempTrend

# como resultado obtenemos 3 capas: spl Trends en ºC/yr;  se Trends standard error; sig Trends statistical significance
# slp simple linear pattern, es decir la pendiente de la regresión linear en cada celda

temp_trend <- tempTrend(yr_sst, th = 10) # 3 nlayers slp, se y significance; 101.640 celdas
temp_trend
## class      : RasterStack 
## dimensions : 180, 360, 64800, 3  (nrow, ncol, ncell, nlayers)
## resolution : 1, 1  (x, y)
## extent     : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
## crs        : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
## names      :     slpTrends,      seTrends,     sigTrends 
## min values : -7.518329e-02,  0.000000e+00,  1.094215e-17 
## max values :    0.07100648,    0.01080910,    1.00000000
plot(temp_trend, main = c("Tendencia temporal de SST (ºC/año)", "SE", "p-valor"))

# otra manera de dibujarlo
plot(temp_trend[[1]], main = "Tendencia temporal de SST (ºC/año)")

plot(temp_trend[[2]], main = "SE Error Estándard")

plot(temp_trend[[3]], main = "p-value")

# Ya solo llegar hasta aquí y obtener la métrica de la tendencia temporal de SST puede sernos muy útil para muchos análisis. 
# esto mismo podríamos calcularlo solo para los meses de primavera, por ejemplo
# Los siguientes pasos son ligéramente más complejos

Cálculo del gradiente espacial

# Calcular el gradiente espacial, promediado para toda la serie temporal. Para más información sobre cómo se calcula, ver la presentación en PowerPoint Métricas climáticas

?spatGrad
spa_grad <- spatGrad(yr_sst, th = 0.0001, projected = FALSE)
spa_grad # obtienes el gradiente y el ángulo, pero solo utilizaremos el gradiente
## class      : RasterStack 
## dimensions : 180, 360, 64800, 2  (nrow, ncol, ncell, nlayers)
## resolution : 1, 1  (x, y)
## extent     : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
## crs        : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
## names      :         Grad,          Ang 
## min values :        1e-04,        0e+00 
## max values :   0.02835585, 359.99913929
plot(spa_grad)

plot(spa_grad$Grad, main = "Gradiente espacial térmico (ºC/km)")

Cálculo de la velocidad climática (gradient-based velocity gVoCC)

# Ahora calcularemos la velocidad climática basada en el gradiente, o también llamado método de las pendientes, porque está basado en el cálculo de las tendencias temporales (tempTrend) a través de las series temporales de temperaturas de cada año, y la pendiente de la recta de regresión en cada celda focal o píxel. 

gvel <- gVoCC(temp_trend , spa_grad)

gvel # dos resultados, la magnitud y el ángulo. Solo utilizaremos la magnitud. 
## class      : RasterStack 
## dimensions : 180, 360, 64800, 2  (nrow, ncol, ncell, nlayers)
## resolution : 1, 1  (x, y)
## extent     : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
## crs        : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
## names      :   voccMag,   voccAng 
## min values : -229.2733,    0.0000 
## max values :  255.9260,  359.9958
gvel_mag <- gvel[[1]] # raster con la magnitud de la gradient-based velocity, voccMag

plot(gvel_mag, main="Velocidad de cambio climático (km/año)") # km/yr

hist(gvel_mag, main="Velocidad de cambio climático")

# otra manera de plotearlo

my.at <- seq(-50, 50, by = 5)
p1 <- rasterVis::levelplot(gvel[[1]], par.settings = BuRdTheme, at=my.at, main = 'VoCC basada en el gradiente (km/año)', margin = FALSE)

gridExtra::grid.arrange(p1)

# comparar con este resultado del creador del paquete: https://github.com/JorGarMol/VoCC/blob/master/images/Fig1.png


# Para guardar y poder abrir en programas como QGIS. Aquí se crearán varios ficheros espaciales en la carpeta de tu directory working
# para saber o recordar dónde se guardarán utiliza lo siguiente:
getwd() # ve al resultado para obt
## [1] "/Users/Marina/Desktop/VICENÇ MUT inicio/OHW Intercoonecta/OHW 2025/Taller intermedio Oct 2025/Marina/Datos espaciales de métricas climáticas"
writeRaster(gvel, "gvel", overwrite=TRUE)
writeRaster(temp_trend, "temp_trend", overwrite=TRUE)
writeRaster(spa_grad, "spa_grad", overwrite=TRUE)

getwd() # ve a esta dirección para obtener los ficheros que has generado
## [1] "/Users/Marina/Desktop/VICENÇ MUT inicio/OHW Intercoonecta/OHW 2025/Taller intermedio Oct 2025/Marina/Datos espaciales de métricas climáticas"

Cortar el póligono para ver solo la zona que nos interesa

# Como el mapa de datos es es muy grande, lo queremos cortar
med_pol <- as(extent(-7.36,7, 30, 42.5), 'SpatialPolygons') # estas son las coordenadas del mediterráneo occidental, latitud, longitud de un polígono. 

# si queremos jugar con las coordenadas, modificamos estos números
# ?extent # más información sobre cómo hacerlo en inglés

crs(med_pol) <- crs(gvel) # para poner el sistema de coordenadas de referencia igual que en el fichero sst_dias, coordinate reference system (CRS) of a raster object.

gvel_med <- crop(gvel, med_pol) # lo corto solo para el mediterráneo occidental, no me interesa el atlántico
gvel_med
## class      : RasterBrick 
## dimensions : 12, 14, 168, 2  (nrow, ncol, ncell, nlayers)
## resolution : 1, 1  (x, y)
## extent     : -7, 7, 30, 42  (xmin, xmax, ymin, ymax)
## crs        : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
## source     : memory
## names      :   voccMag,   voccAng 
## min values :  1.295291,  3.607703 
## max values :  51.35228, 359.43534
plot(gvel_med[[1]])# ?crop del paquete raster

# puedes usar Ctrl+Shift+C para agregar o eliminar comentarios en varias líneas de código.
# puedes usar el mismo atajo de teclado para eliminar el indicador de documentación #